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We present a microscopic calculation of the Gilbert damping constant for the magnetization of 
a two-dimensional spin-polarized electron liquid in the presence of intrinsic spin-orbit interaction. 
First we show that the Gilbert constant can be expressed in terms of the auto-correlation function 
of the spin-orbit induced torque. Then we specialize to the case of the Rashba spin-orbit interaction 
and we show that the Gilbert constant in this model is related to the spin-channel conductivity. 
This allows us to study the Gilbert damping constant in different physical regimes, characterized by 
different orderings of the relevant energy scales - spin-orbit coupling, Zeeman coupling, momentum 
relaxation rate, spin-momentum relaxation rate, spin precession frequency - and to discuss its 
behavior in various limits. Particular attention is paid to electron-electron interaction effects, which 
enter the spin conductivity and hence the Gilbert damping constant via the spin Coulomb drag 
coefficient. 
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I. INTRODUCTION 

The Gilbert constant characterizing the damping of 
magnetization precession is one of the important phe- 
nomenological parameters that describe the collective 
magnetization dynamics of ferromagnets J^i 3 -. It is 
an essential input of the Landau-Lifshitz-Gilbcrt equa- 
tion (LLG) for magnetization dynamics 3 - and as such 
is widely used in the analysis of magnetization rever- 
sal processes, which are crucial to magnetic recording 
technologies 4 .. Recently, hybrid systems of fcrromagnets 
and normal metals have also attracted considerable at- 
tention, in the context of the enhancement of Gilbert 
damping at the ferromagnet-normal metal interfaces'* 5 ^ 
Despite tremendous efforts to elucidate the nature of 
Gilbert damping in bulk fcrromagnets, however, the mi- 
croscopic processes responsible for the observed ferro- 
magnetic relaxation in real materials are still not fully 
understood. The matter is complicated by the po- 
tential relevance of a number of different mechanisms 
involving eddy currents, magnetoclastic coupling, two- 
magnon and nonlinear multimagnon processes, extrin- 
sic and intrinsic spin-orbit (SO) coupling of itinerant 
electrons^^°^i 12 i 13 i 14 i 15 i 16 t 17 . 18 i 19 i 20 Without delving 
into a detailed discussion of various points of view, we 
note that a certain level of consensus has been reached 
in the recent theoretical literature about the central im- 
portance of SO interaction of some form in conducting 
fcrromagnets, although with a limited and rather indi- 
rect experimental support at present i 21 ' 22 

Even restricting the attention to the SO-based mecha- 
nisms of Gilbert damping, however, the myriad of rel- 
evant energy scales (namely, ferromagnetic resonance 
frequency, ferromagnetic exchange energy, intrinsic SO 
splitting, impurity scattering rate, and spin dephasing 
due to magnetic or SO disorder) has led different au- 
thors to make qualitatively different predictions for the 



Gilbert damping, for example, with regard to its depen- 
dence on the disorder strength ] 7 ! 13 : 15 ' 16 : 20 ' 23 . Further- 
more, electron-electron (e — e) interactions have mainly 
been discussed in the mean- field spirit, as the source of 
the exchange field in itinerant-electron ferromagnets. 

In this paper, we set off to formulate a microscopic the- 
ory of the Gilbert damping constant, which wc then apply 
to a simple model, where the competition of different en- 
ergy scales can be studied, comparing various points of 
view and also generalizing and extending the results ex- 
isting in the literature. We wish to consider the interplay 
of magnetic spin splitting, intrinsic SO strength, disor- 
der scattering, as well as the strength of e — e interactions 
in a single self-contained model without introducing any 
phenomcnological dephasing or relaxation parameters. 
In particular, we find an intricate relation between the 
spin-drag correlations induced by e — e interactions and 
the Gilbert damping in two-dimensional electron liquids. 
We hope our discussion will bring the community a step 
closer to understanding the relevant microscopic mech- 
anisms responsible for the Gilbert damping observed in 
conducting ferromagnets and spin-polarized systems. 

SO interactions, especially in narrow-band semicon- 
ductors, have recently received a great deal of attention 
in the context of spin-based electronics (i.e., spintronics). 
In particular, theoretical proposals to manipulate spins 
by means of intrinsic SO interactions, without the use of 
magnetic fields or magnetic materials, especially in the 
so-called spin Hall configuration ) 24 ' 25 have unleashed a 
wave of theoretical research as well as experimental ef- 
forts to measure the effectj 26 ' 27 : 28 It is worthwhile not- 
ing that these activities have shown the need to better 
understand the fundamental aspects of intrinsic SO cou- 
pling and its experimental manifestations. The role of 
e — e interactions, furthermore, remains a relatively un- 
explored territory with regard to its interplay with in- 
teresting topological properties brought about by the SO 



coupling, which lie beyond the conventional theories of 
electron liquids. In addition, the recent discovery of 
the Aharonov-Casher phased in narrow-band semicon- 
ductors shows that a study of the combined response to 
magnetic and SO fields in a regime of weak and strong SO 
interactions is of urgent importance. It is closely related 
to the topic considered here, i.e., the spin response and 
relaxation in magnetized system with SO interactions. 

In this paper, we study an interacting electron liquid 
magnetized by an external magnetic field in the presence 
of SO interactions. First, by comparing the macroscopic 
Landau-Lifshitz-Gilbcrt (LLG) theory with a microscop- 
ically derived expression for the transverse spin suscep- 
tibility, we find the following relationship between the 
Gilbert damping constant and the torque-torque correla- 
tor (see Eq. 32 below): 
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where g is the gyromagnetic ratio, A is the dimensionless 
Gilbert coefficient, Mq is the equilibrium magnetization, 
V the volume of the system, and f the torque induced 
by SO interactions. To derive Eq. |T]) we must neglect 
certain contributions of order higher than the fourth in 
the strength of the spin-orbit coupling: we will argue that 
this procedure is justified, not only at zero frequency but 
also at frequencies close to the ferromagnetic resonance. 

Next we focus on the specific case of a two-dimensional 
electron gas (2DEG) with Rashba SO interactions in- 
duced by a uniform electric field perpendicular to the 
plane of the electron gas. In this case, we show that 
the Gilbert damping constant can be expressed in terms 
of the spin-channel conductivity: the same result holds 
in a three-dimensional electron gas, provided the electric 
field responsible for the spin-orbit interaction is uniform 
in space. For the isotropic case, i.e., when the magnetic 
field, the magnetization, and the SO-inducing electric 
field are all perpendicular to the direction of the 2DEG, 
the result is 



A = 



a 2 E F T^ea s \\{{)) 
hp a D 



(2) 



where Ep is the Fermi energy, r is the momentum relax- 
ation time due to electron-impurity scattering, a is the 
Rashba SO-coupling strength normalized by the Fermi 
velocity Vp, p is the degree of spin polarization, op, is 
the Drude conductivity, and c s ||(0) is the longitudinal 
spin-channel conductivity. 

The problem now shifts to the calculation of the spin- 
channel conductivity. We show that this can be done 
exactly (subject to the usual weak-disorder assumption 
EpT ^> 1) in the isotropic non-interacting case, with lad- 
der vertex corrections playing an essential role in ensuring 
the correct behavior in the limit of strong SO coupling, 
e — e interactions can then be introduced using the for- 
malism described in Refs. I30ll3ll. The final result for the 



isotropic case at zero frequency has the form: 

x _ pa 2 E F T 1 + i(fir) 2 (l + cos 2 <5) l+p 2 1 r 
~ % 2 +a 2 )cos 2 <5+±(fV) 2 (l + cos 2 c5) 2 1+ 7 T ' 

(3) 

where Q is the spin splitting of electronic states due to 
the combined action of the external magnetic field and 
the SO effective magnetic field, cos 8 is the ratio of the ex- 
ternal magnetic field to the total effective field fl, and 7 
is the spin-drag coefficient. This expresses the Gilbert 
damping function of disorder strength, magnetic 

field, SO and e — e interactions. 

One of the gross features of Eq. is a rough scaling of 
A with t, which is particularly evident in the clean limit 
t — ► 00 . We will see, however, that Eq. §5§ describes 
in general a nonmonotonic dependence of A on the scat- 
tering rate, due to the interplay between 1/t and other 
energy scales. For example, in the regime of strong e — e 
interactions and weak polarizations, 1 < jt < l/p 2 , the 
Gilbert damping A becomes independent of r and scales 
as I/7. 

It is important to note that Eq. ([3]) holds in the low- 
frequency limit, uj — > 0, which means in practice uj much 
smaller than all the other energy scales. For this rea- 
son the divergence of A as 1/t — > in Eq. ([3]), which 
was obtained after first taking the uj — > limit, is not 
physically consequential: taking the 1/t — > limit with 
a fixed uj results in a vanishing damping as it should^. 
This is clearly seen by extending the calculation of A to 
finite frequency, which can be done with little extra ef- 
fort in the uj <C Ep regime. This important extension of 
Eq. is presented in full in Eq. 59 below. Among other 
modifications we get a factor 1+ ^, T )3 to multiply Eq. (J3)), 
ensuring the correct behavior of A for r — > 00 at finite uj. 

The rest of the paper is organized as follows. In Sec. ITU 
we discuss the LLG phenomenology, which is compared 
to the microscopic calculation in Sec. IIIH allowing us to 
express the Gilbert damping constant in terms of the 
torque-torque correlator. In Sec. IIV Al we present the 
calculation of the Gilbert damping constant in the pres- 
ence of the Rashba SO and e — e interactions for an 
isotropic case, i.e., when the magnetic field is perpendic- 
ular to the 2DEG plane. In Sec. IIV Bl the same model is 
treated for the anisotropic case, where the induced mag- 
netization is tilted away from the direction perpendicular 
to the 2DEG. Our summary, conclusions and specula- 
tions are presented in Sec. [V] Technical details of the 
calculations are supplied in the Appendixes. 



II. LLG EQUATION 

Studies of the magnetization dynamics in ferromag- 
netic materials often start from a phcnomcnological equa- 
tion of motion for the magnetization 



M{r,t) =-<7(s(r,t)>, 



(4) 
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where (s(r, t)) is the expectation value of the spin-density 
operator at position r and g is the gyromagnetic ratio: 
for free electrons in vacuum, g — e/mc where e is the 
absolute value of the electron charge. An important as- 
sumption in the LLG phenomenology is that the magni- 
tude of the magnetization is fixed at a constant value Mo 
(determined by the minimization of the free energy), i.e., 
we have 



from the simplest model, and then try to generalize our 
conclusions to more realistic situations. Note that al- 
though in our model we consider SO coupling stemming 
from the relativistic corrections to the Schrodinger equa- 
tion, a kind of SO coupling would also be possible in the 
presence of an inhomogcncous crystal magnetic field. 



M = M e M , 



(5) 



where eM is the unit vector specifying the direction of 
the magnetization. The physical reasoning^ underlying 
this picture assumes a very large (essentially infinite) 
longitudinal spin stiffness which prevents changes in the 
magnitude of M, while one focuses on the softer trans- 
verse modes - essentially ferromagnetic Goldstone modes 
with a finite gap in the presence of intrinsic or extrinsic 
anisotropics. The direction of the magnetization §m thus 
provides the relevant degrees of freedom, whose time evo- 
lution is determined by the LLG equation 



= .gHoff x e M + A - eM x 



dt 



dt 



(6) 



where H e g is the effective magnetic field defined as the 
functional derivative of the equilibrium free-energy den- 
sity /(M), regarded as a functional of the magnetization 
M, with respect to its argument: 



Heff(r,i) = - 



*/(M) 
£M(r,t) 



(7) 



A is a symmetric 2x2 matrix acting in the two- 
dimensional space perpendicular to the magnetization, 
and the dot denotes the usual matrix product [note that 

the possible antisymmetric component of A would only 
lead to a renormalization of the gyromagnetic ratio in 
Eq. 0.] The first term on the right-hand side (r.h.s.) 
of this equation produces a coherent processional motion 
of the magnetization, conserving the total free energy, 
while the second term, known as the Gilbert damping, 
is responsible for the relaxation toward equilibrium, i.e., 
the lowest free-energy state. Both terms vanish at equi- 
librium. Notice that the LLG equation is not restricted 
to the linear-response regime, as e^j is allowed to wander 
away arbitrarily from the equilibrium orientation. 

The primary focus of this paper is on the microscopic 
calculation of the Gilbert damping matrix A • We con- 
centrate exclusively on homogeneous systems, that is the 
situation in which M is independent of position and de- 
pends only on time. In such a case, it is immediately 
obvious that any damping must arise from the coupling 
between the spin and the orbital degrees of freedom. This 
is so since the nonrelativistic e — e interaction is invari- 
ant under a uniform rotation of all the spins, and thus 
cannot be involved in the damping of such a rotation. It 
will therefore be essential to include the SO coupling in 
our microscopic calculation, and we will do this starting 
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FIG. 1: The directions of the magnetic and Rashba fields in 
relation to the 2DEG plane, the model which we discuss in 
Section IV. Also, visualization of directions of equilibrium and 
induced magnetization presented in Eqs. 8-14. 



We can calculate A by applying Eq. ((BJ) to a linear- 
response problem and comparing the resulting suscepti- 
bility with the susceptibility obtained from the Kubo for- 
mula. Our linear-response problem is depicted in Fig. [1] 
We perturb the system with a small magnetic field Hi 
perpendicular to the direction of the equilibrium magne- 
tization, which we denote by §3 to distinguish it from the 
actual direction of the magnetization ejj. Hi is uniform 
in space and periodic in time, with an angular frequency 
uj (this field must not be confused with the static field 
Ho which may be present at equilibrium). In response 
to Hi, the direction of the magnetization becomes time 
dependent and will be written as 



e M (i) = e 3 + Se M (t) 



(8) 



where <5eM is a small vector perpendicular to §3, which 
therefore does not change the normalization of &m to 
first order. Adding the external field Hi to H e g and 
substituting Eq. ||5J) into Eq. ([6]), we obtain the equation 
of motion for <5eM to first order in Hi: 



dSe M 



gr(Hi + H cff ,i) x e 3 + A • e 3 x 



<9<5eM 



(9) 
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To the same order, the effective field is 

_J_ d 2 f(M) 
M§ de M de M 



H 



cff.l 



M Se M 



= - S -Mo5e M , (10) 
where we have defined the transverse spin stiffness matrix 



j_ d 2 f(M) 

Mq <9e M <9e M 



(11) 



a 2 x 2 matrix in the plane perpendicular to 63. 

Upon taking a Fourier transform with respect to time, 
the linearized LLG equation can be rewritten as follows: 

-iu)(l + 7 ■ A) + gM 7 ■ S -6e M = 5 ~ • H x , (12) 

where 1 is the 2x2 identity matrix and 

1 



-1 



(13) 



e = — e . The transverse spin susceptibility X con- 
nects the transverse magnetization Mi = AlgSe-M to the 
external field Hi: 

Mi(w) =X (w)-Hi(w). (14) 
Comparing this with Eq. (112p . we easily obtain 



X H= 77 +S, 

gM 



(15) 



which is the most general form of the anisotropic trans- 
verse response, to the first order in frequency. Notice 

that X (0) =S, showing that the transverse stiffness 
is indeed the inverse of the static spin susceptibility. In 

an isotropic system, S vanishes (the free energy does 
not depend on the direction of the magnetization). If a 
static external field Ho = HQ&3 is applied to such a sys- 
tem, then the stiffness matrix becomes S= {Hq/Mq) 1 . 
Anisotropics due to SO interactions will show up as ad- 
ditional contributions to S- 

Needless to say, the theory described in this section is 
valid only at low frequencies (in comparison to relevant 
microscopic energy scales), when the frequency expan- 
sion of the inverse susceptibility, such as Eq. (fl"5"j) . can be 
truncated at the linear term. 



III. MICROSCOPIC THEORY 

A microscopic expression for the homogeneous trans- 
verse spin susceptibility is given by the Kubo formula 



XyM = — 



dt([Si{t),sM\) 



, (16) 



where Si(t) is the Hcisenberg operator of the total spin 
of the system, the angular bracket denotes the equilib- 
rium average, V is the volume of the system, and i, j 
are Cartesian indices in the plane perpendicular to the 
magnetization (they can take up the values 1 or 2). The 
lengthy expression inside the square brackets in the above 
equation (i.e., the retarded spin response function) will 
be abbreviated from now on as —((Si'jSj}}^, so we have 



Xij(«>) = -y((Si;Sj)) 



(17) 



(The magnetization is related to the average spin by M = 
—g(Si)/V.) Notice that this expression is dimcnsionlcss 
in three dimensions and has the units of length in two 
spatial dimensions (in the cgs units). 

In view of the fact that Gilbert damping arises ulti- 
mately from torques on the individual spins, it is conve- 
nient to rewrite Eq. (|17j) in terms of the time derivatives 
of S. The time derivative of the spin is conveniently 
written as the sum of two parts: a free precession term 
at frequency loq about the magnetization axis §3 and a 
residual torque r due to the SO interactions: 



S = 0^063 x S 



(18) 



Here, uj^ = gH' , where H' Q is the projection of the static 
magnetic field along §3, while f is a model-dependent 
torque, which is proportional to the spin-orbit coupling 
constant. 

It is very convenient at this point to rewrite Eq. (|17p 
in terms of the ((t;t)) correlation function. The fact 
that t explicitly contains the spin-orbit coupling con- 
stant a will allow us to calculate the correlator, to order 
a 2 , without including the spin-orbit interaction in the 
hamiltonian: this will be a major simplification. Fur- 
thermore, this transformation will enable us to connect 
the Gilbert constant to the spin conductivity, which is 
particularly helpful for the inclusion of e — e interactions. 

To express the spin susceptibility in terms of the 
((t; t)) correlator we make use of the identity 



{{A;B)) U = 



((A;B))„ | ([A,B]) 
(A , ([AB]) 



where 



A = -[A,H] 
in 



(19) 



(20) 



is the time derivative of A and H is the unperturbed 
Hamiltonian. 

Combining Eqs. (Til?)) and p8[) . we easily obtain 



-i U ((S i ;S j )) u = ((S i ;S j )) u - l {[S i ,S j ]) 
eij{Sa) -v eik((Sk;Sj})u> + {{Ti;Sj})<j, (21) 
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where a sum over the repeated index (k in this case) is 
understood. Solving for ((Si] Sj))^, we get 



and 



(22) 

In the absence of the SO coupling, the torque r vanishes 
and the second term on the r.h.s. of Eq. (f2"2"P is absent. 
Then the susceptibility is easily found to have the form 



X$ ) (u)=gM '- 



Notice that its inverse 



iuj e +ujq 1 



gM 



(23) 



(24) 



is in a perfect agreement with Eq. I|15p. if in the latter 

we set A = and S= (uo/gMo) 1 . 

Upon applying again Eq. (TH)]) to (("?&; we get 

iuiSij + u> eij 



-([f fe) 5 J ]) + ((f fe ;f i )) a 



(25) 

Putting this in Eq. (|22p and making use of Eqs. (jTTJ) 
and (|23p . we arrive after some algebraic manipulations at 
the following compact formula for the spin susceptibility 
matrix: 



~ «(0) «(0) e • r M- e ~(°) , s. 
X (o>) =X (w)+ X (u3) —5-7 X (w) 



M 2 U 



where the matrix r (w) is defined as follows: 



(26) 



(27) 



The main dynamical quantity on the r.h.s. here is the 
torque-torque response function ((Ti; Tj})^. The rest is 
ground-state averages that we will not be interested in. 
We can now express the Gilbert damping and stiffness 

matrices, A and S, which appear in the phenomenologi- 
cal equation (|15p , in terms of microscopic response func- 
tions. First of all, consider the formal limit of weak SO 
interactions. In this limit, r, which is proportional to the 
square of the SO coupling constant, is assumed to give 
a small correction to the susceptibility (we will comment 
below on the conditions of applicability of this regime). 

Upon inverting Eq. (|26[) to first order in r , we obtain 



x M ~ [x 



r M- e 



(28) 



Comparing this with Eq. (|15p . and taking into account 
Eq. (HID, we get 
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r («)■ 



MoV 



(30) 



where the subscript sym denotes the symmetric part 
of a matrix. In the first of these two equations, S so 
is the SO contribution to the stiffness matrix, i.e., 

Sso=S -~{Hq/M ) 1, which reflects SO induced mag- 
netic anisotropy. The second equation is the main result 
at this point. It expresses the Gilbert damping matrix 
in terms of the zero-frequency slope of the spectrum of 
the torque-torque response function. Notice that this 
spectrum has no contribution from the second term on 
the r.h.s. of Eq. (|27p . which is purely real (the commu- 
tator of two hermitian operators is anti-hcrmitian) . In 

general, the imaginary part of e ■ r (w)- e will have 
both symmetric and antisymmetric components: the lat- 
ter has been interpreted as a Berry curvature correction^ 2 - 
to the adiabatic spin dynamics [also, see the note below 
Eq. ([7]), where we point out that inclusion of the antisym- 
metric component would lead to a renormalization of the 
g-factor]. The symmetric component is purely diagonal 
in the limit of weak spin-orbit coupling because in this 
limit the system is isotropic under rotations about the 3 
axis. (In general, one can still diagonalize the symmetric 
component of the spectrum by a suitable choice of the 
axes in the 1 — 2 plane - a transformation that does not 
affect the antisymmetric component. We will not con- 
sider this complication here) . Thus we conclude that the 
Gilbert damping matrix is a diagonal matrix of the form 



where 



and 



Ai 
A 2 



1 9 v ex. ((h;r 2 ))i, 

M = — t~t~t7 hm xsm 

M V w-^o lu 



M V ^0 



(31) 



(32) 



(33) 



The above formulas have been derived from a first or- 
der expansion in r, which is justified for u <C ujo when 
the spin-orbit interaction is weak. At first sight, how- 
ever, the approximation seems to break down completely 

when uj ~ ujq (the interesting region of the ferromagnetic 

«(0) 

resonance) because X has a pole at u> = ujq. But this 
conclusion is too hasty. To show this, we first notice that 
the exact inversion of X from Eq. ([26]) gives 



S so 



r(o> 



M 2 U 



(29) 



,(0) 



x H = [x ]~V)- s M 



(34) 
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where the "self-energy" £ (cu) is given by 



£ (w) = 



« <->(0) 

1 + X (w) 



Af 2 V 



,(0) 



(35) 

For cj < cj , X (w) reduces to Sq 1 = a ^ s -, which 
is just the inverse of the stiffness to zero-th orcler in the 
spin-orbit interaction. Then, from the above formulas 
and making use of Eq. ([29]) we see that the zero- frequency 

Gilbert matrix (proportional to 3to £ (uS)/uS) is given by 
Eq. (|30p times a renormalization factor 



(36) 



We will estimate this renormalization below for a specific 
model, and find it to be very small. 

An important observation is that, unlike the function 

<-(o) 

X (cl>), the self-energy is well-behaved at to = uq. In- 
deed from Eq. (|26[) we see that T M must vanish at oj 

in order to cancel the pole of X and replace it with a 
roughly Lorcntzian peak at the true ferromagnetic reso- 
nance. Furthermore, the rate at which r (to) vanishes 
for lu tending to luq must be such that the denomina- 



tor of Eq. ([35]). i.e., 



f 



.(0) 

x M 



■rM- 



must also 



vanish at ujq - otherwise the pole at luq would survive. 
(Note: we say that an operator vanishes or has a pole in 
the sense that its projection along the eigenvector of the 
ferromagnetic resonance, (l,i) or (1, — i), has a zero or a 

pole). Therefore £ (u>o) is finite and amenable to a per- 
turbativc treatment in which we retain only the leading 



term £ (ui) 
formula 



■rM-' 



This leads to the more general 



A M 
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M V 



sin- 



rM- 



(37) 



which is applicable at finite frequency, as well as zero 
frequency. Furthermore, we will show below that this 

approximation for £ generates results that are in perfect 
agreement, in the non-interacting case, with the direct 
diagrammatic calculation of the spin susceptibility, even 
when a is not small. 



a linear SO interaction induced by an electric field along 
the z axis perpendicular to the plane of the electrons, 
and a magnetic field H = where the unit vector §3 

forms an angle with the z axis and lies in the z — x 
plane. The hamiltonian for the model is 



n 



^- + V(r n ) + uj S n • e 3 - ^(e z xp)-S, 
Ira a 



2 <-f e 6 |r„ - iv I 



(38) 



where p„, f„, and S n are, respectively, the momentum, 
the position, and the spin operators of the n-th electron, 
V(r) is the random electron-impurity potential, a is the 
Rashba velocity, which controls the strength of the SO 
coupling (a is proportional to the electric field in the z 
direction) , m is the effective mass of the electrons and ef, 
is the background dielectric constant. The explicit form 
for a in a typical direct gap semiconductor (say GaAs) is 



a = vp(eEaB) 



2P' 2 
3m 



1 



1 



(E g + A so y 



V2r s 



(39) 

where vf is the Fermi velocity, E is the magnitude of the 
electric field in the z direction, as is the Bohr radius, P is 
the matrix element of the momentum operator between 
the conduction and the valence band at the zone center, 
Eg is the gap, A so is the SO splitting of the lowest hole 
band, and r s is the average distance between the electrons 
in units of as- 

Notice that we have already assumed that the direction 
of the external magnetic field coincides with the direc- 
tion of the equilibrium magnetization (§3 , by definition) . 
This is not generally true in the presence of SO coupling 
(except for special cases such as O = 0), but the an- 
gle between Ho and Mo is of order ( h E k J \ and will 

be neglected henceforth. The magnitude of the magneti- 
zation, neglecting Coulomb interactions, is simply given 

by M = (sf) N(0)H where N(0) = ^ is the free- 
particle density of states in two dimensions. 

The SO torque, defined in Eq. (fT8|) is straightforwardly 
calculated to be 



(40) 



IV. SPIN RESPONSE OF A RASHBA 2DEG 

In this section we introduce a simple model that allows 
us to calculate the torque-torque response function and 
hence the Gilbert damping and the static stiffness from 
first principles. The model is illustrated in Fig. 1 and 
consists of a two-dimensional electron gas (2DEG) with 



In order to facilitate the calculations in the limit of weak 
SO coupling it is convenient to express the torque in 
terms of the components p x and p y of the momentum 
in the plane, as well as the components Si, S2 and S3 of 
the spin in the coordinate system shown in Fig. 1. This 
is advantageous when the only significant anisotropy in 
spin space is caused by the magnetic field in the §3 di- 
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rection. The result of this rewriting is 

fi = ah^ 1 J2 n (s n3 p nx - S n2 p ny sine) 

T 2 = ah" 1 J2 n (Sn,3pny COS + SnlPny sin 9 ) (41) 
f 3 = ah' 1 J2„ \ ~SnlPnx - S n2 pnyCOsQ) 



Notice that these equations do not involve p z - not just 
because our model system is two-dimensional but, more 
fundamentally, because a motion along the direction of 
the electric field does not produce spin-orbit coupling. 

The commutators [r,,^-] that appear in Eq. (|27[) can 
also be straightforwardly calculated: 



i[fi, Si] 
-i[T2, S 2 ] 
-i[h, S 2 ] 
-i[f 2 ,Si] 



a ( S n sp n y sin + 52 nPr, 

n 

-OL [SnlPny COS - S n SPny Sm 0J 
n 

nlPnx 

n 

Ot X! Sn2Pny COS © . (42) 



Their expectation values in the ground state are of sec- 
ond order in a and are closely related the change in en- 
ergy of the electron gas due to the SO interaction (i.e. 
they control the magnetic anisotropy induced by the SO 
interaction) . 

The fact that in this model the torque is linear in 
both spin and momentum, i.e. proportional to the spin- 
current, allows us to establish an exact connection be- 
tween Gilbert damping and spin-channel conductivity, 
defined in terms of spin-current-spin-current response 
functions. Unfortunately, this relation does not hold in 
more realistic models of the SO coupling, in which the 
electric field that is responsible for the SO coupling de- 
pends on position rather than being constant. But we 
will argue in the concluding section that a similar rela- 
tion, involving a spatially varying SO coupling constant 
and a local spin conductivity can be justified if certain 
conditions arc met. 



A. Isotropic case 

1. Non-interacting 2DEG 

Let us begin with the simplest case in which the mag- 
netic field and the magnetization are parallel to the z 
axis (0 = 0). The system is invariant for rotations about 
the z axis, which coincides with the 3 axis, so we have 



({n;fi))a, = ((T 2 ;f 2 )) G 



a 



((^2 SnSPnx; ^ SnZPnx))u ■ (43) 



Let us introduce the longitudinal spin-channel conduc- 
tivity, a s \\ , as the constant of proportionality between the 



spin-current jf — jj_ and an electric field E s , which acts 
with opposite signs on the spin up and the spin down 
components of the electron liquid, 35 f and J, denote the 
directions parallel and antiparallel to the magnetic field 
(the §3 axis). <7 s ii is related to the spin-current spin- 
current response function by the Kubo formula 

n Sn3Pnx 



mrV 



h 2 o 



3?e<T s || (uj) 

. (44) 

We now discern a simple relationship between the imagi- 
nary part of the torque-torque response function and the 
spin-channel conductivity, namely 



1 3m T (w) 



4e 2 



-$tea s ii(u)) 1 



(45) 



Substituting this in Eqs. ([32H33]) we arrive at the fol- 
lowing expression for the Gilbert damping constant: 



A = Ai 



A, 



9 2 

m ga 



4A/ r 



Ke<M(0) . 



(46) 



Let us emphasize that the relation between equilibrium 
magnetization and external field in this model is not af- 
fected by SO interactions as long as both spin-orbit split 
bands are occupied. The e-e interaction, on the other 
hand, would reduce the value of u>q required to produce 
a given Mo, since part of the magnetization would then 
arise from the exchange interaction. However what will 
be needed in the subsequent calculations is not loq per se, 
but the Zecman splitting that it produces on single par- 
ticle energy levels. The latter is brought back to the non- 
interacting value once exchange is included in the single 
particle energy levels. So it seems permissible to ignore 
the e-e correction to uq, as long as we do not include 
e-e effects in the spectrum of single particle excitations. 
Finally, it can be shown that the spin-orbit anisotropy 

stiffness S so for this model begins with terms of order 
(q;/uf) 4 and is therefore negligible for typical values of a 
and vp. 

Having considered all this we can express A in the ele- 
gant form 



A 



a 2 E F r Ketr a ||(0) 



hp 



(TD 



(47) 



where a is the Rashba velocity in units of the Fermi ve- 
locity 



<7d is the usual Drude conductivity 



id 



(48) 



(49) 



r is the transport scattering time, Ep is the Fermi energy, 
and p is the degree of polarization of the electron gas 



V 



hu>o _ 2M 
2Ef ghn 



(50) 
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The next task is the calculation of the spin-channel 
conductivity. In the absence of Coulomb interaction this 
can be done even without assuming that a is small com- 
pared to p and H/EpT. We supply the details of the 
calculation in Appendix A. The result is 



» Krt 2e l + i(^) 2 (l+COS 2 5) 

iRea s « (0) = <7d cos o- 



cos 2 <5 + i(fir) 2 (l +cos 2 8f 



where we have defined 



and 



si = 



, w 
cos = — = 



a 2 k F 



and correspondingly: 



Vi 2 



(51) 



(52) 



(53) 



(54) 



We have assumed, for simplicity, that the transport scat- 
tering time t is the same for up- and down-spin electrons. 
Notice that, in the absence of SO coupling, one recovers 
<7 s |i (0) = <7d as expected since up- and down-spin com- 
ponents are completely decoupled and have the same mo- 
bility. If, on the other hand, the external field frequency 
luq is set to zero when a is still finite we get c s ||(0) = 0. 
Again, this is not surprising in view of the fact that the 
SO interaction causes a steady precession of the spin in 
a plane perpendicular to p, effectively suppressing the 
average z-component of the spin. Putting Eq. (fSTj) back 
into Eq. (|4T[) we arrive at 



E F Tsin 2 5 l + ±(Slr) 2 (l + cos 2 <5) 
' H cos 2 <5 + ±(Slr) 2 (l + cos 2 5f 



(55) 



where the subscript ni stands for non-interacting. Notice 
that this expression vanishes in the limit of zero magne- 
tization, i.e. for p — > 0. 

The Gilbert damping constant, as given by Eq. (|55p is 
plotted in Fig. 2 as a function of two parameters: Sir - 
measuring the effectiveness of spin precession during an 
elastic mean free path r - and akp/uj = a/p the ratio of 
the SO effective magnetic field to the external magnetic 
field. Clearly Sir 3> 1 corresponds to the "clean" limit 
and Sir < 1 to the "dirty" limit. 

Fig. 3 shows the two-dimensional cross-sections of 
Fig. 2. In Fig. 3a we see the Gilbert constant plotted as 
a function of the spin-orbit to polarization ratio for two 
values of Sir, one in the clean limit and one in the dirty 
limit. In the weak SO coupling regime (akp/tuo -C 1) the 
Gilbert constant grows as (afc^/wo) 2 , while for strong SO 
coupling it saturates and begins to decrease linearly be- 
yond a certain value akp/ujQ. The quadratic increase of 
A for small a is easily accounted for by the growth of the 
sin 2 S factor in Eq. ([53)) - sec Eq. (f5"4"|) . For large a, on 




FIG. 2: Gilbert damping constant as a function of SO to mag- 
netic field ratio akp/uJo and disorder scale fir. We assumed 



U2D 



10 cm and the polarization p=0.1 



the other hand, sin<5 approaches the maximum value of 
1 and A becomes essentially proportional to r. The de- 
crease in A in this regime reflects the decrease in r along 
a curve on which SI grows while Sir remains constant. 

Fig. 3b shows the Gilbert damping as a function of Sir 
for two different values of akp/uQ. For weak SO coupling 
(lower solid line) the Gilbert damping constant is propor- 
tional to t in both the clean and dirty limits. This is not 
surprising given that in this regime the Gilbert damp- 
ing is essentially proportional to the Drude conductivity. 
For strong SO coupling, however, we observe an interest- 
ing non-monotonic behavior of A (upper solid curve) in 
the transition region between the dirty regime (Sir <C 1) 
and the clean one Sir ^> 1. This region is defined by 
the inequality 2 cos 5 < fir < \/2, and obviously appears 
only if cos (5 is sufficiently small, i.e., if the SO coupling 
is sufficiently strong. 

There is a limit, however, on how large the spin-orbit 
coupling can be made for a given value of Sir. Recall that 
our treatment of disorder is justified for EpT ^> 1 and 
becomes uncontrolled when this inequality is violated. 



Since Eft 



we see that S!r/2p must be 



2p^l + (afc F /w ) 2 

larger than or, equivalently, Sir > 2a if EpT ^ 1 is 
to be satisfied. In Fig. 3b the regions of validity of this 
condition, for three different values of a, extend to the 
right of the three vertical dashed lines. We see that the 
non-monotonic behavior for larger a occurs well within 
the region in which our treatment of disorder is justified 
and can be observed if a is much larger than p, so that 
cos S is small. 

In a naive microscopic picture the Gilbert damping 
could be obtained by substituting the D'yakonov-Perel 
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We note, however, that the full structure of Eq. ([55)) can- 
not be completely captured by the D'yakonov-Perel spin 
relaxation rate (f56]) . In particular in the clean limit the 
combined Eqs. ([55)) - ([57)) would lead to the scaling of 
A as 1/t, which tends to zero, while our Eq. ([55)) gives 
scaling of A as r, which tends to infinity. 



\<x/u F = 0.1 fix =25 



\a /u F = 0.01 
fit =25 



0.0 0.2 0.4 0.6 0.8 1.0 

polarization, p 



FIG. 3: (a) The Gilbert damping constant as a function of 
SO to magnetic field ratio for two different values of disorder 
scale fir corresponding to the clean and dirty limits. The 
magnetic field and as a consequence the polarization is fixed, 
while the strength of spin-orbit interactions and r changes so 
Sir = const. The curves for p — 0.1 and p = 0.0005 and the 
ratio Q,t = 0.1 are indistinguishable, (b) The Gilbert damp- 
ing as a function of disorder scale for two different values of 
afe^/wo corresponding to the strong and weak SO couplings. 
The magnetic field and so the polarization and the spin orbit 
coupling strength are fixed. The condition of validity of our 
calculations Eft 3> 1 is on the right to the dashed line. The 
curves for p — 0.1 and p = 0.0005 and the ratio akp/uJo = 5 
are indistinguishable. 



spin relaxation rati 



,33 



FIG. 4: The Gilbert damping constant as a function of po- 
larization for two different values of renormalized spin-orbit 
coupling a = ol/vf- 



Fig. 4 shows the Gilbert damping as a function of po- 
larization for two different values of spin-orbit coupling. 
One can see that Gilbert damping increases with p for 
p < a and decreases with further increase of p above 
a. Tunable magnetization damping is very desirable 
since it allows to reduce post-switching magnetization 
precession^. Fig. 4 shows the possibility of tuning of the 
Gilbert damping by changing temperature in a polarized 
2DEG with spin-orbit interactions. 



1 



2r{aki 



1 + (2rafc F ) 2 



(56) 



into the phcnomenological Bloch equation for spin dy- 
namics. This results in the formula 



The uj dependent spin-channel conductivity, can be ob- 
tained through the replacement of ^ by - — iu (see also 
appendix A) and has a form: 



A 



1 + {uj T s y 



(57) 



which reproduces Eq. ()55[) either in the weak-SO limit 
with akp <C 1/t,u>q or in the absence of the applied 
field, Wo = 0, where Eq. ([5"6")) becomes the correct Bloch 
spin relaxation rate of the unpolarized Rashba 2DEG^. 



!RecT s |[ (ui) = <td 



1 + (wr) 2 



1 + (^r) 2 [isin^+ (cosS + u/tt) 2 ] 
cos 2 S + (Or) 2 [i sin 2 S + cos <5(cos 6 + 



(58) 
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which gives 



Ani(w) = P 



a 2 E F r 



h(p 2 + a 2 ) 1 + (wt) : 

\2rl 



1 + (fir) 2 [± sin 2 S + (cos <5 + w/ft) 2 ] 
cos 2 (5 + (fir) 2 [i sin 2 <5 + cos <S(cos S + uj/U)} 



(59) 



These equations reduce to Eqs. ([51]) and (|55|) respectively 
if we set a; = 0. 

Notice the very different limiting behavior of A as a 
function of the disorder strength for lo = and u> finite. 
In the first case, A tends to infinity as the scattering time 
t increases (i.e. for decreasing disorder). At finite u> 
however, the Gilbert constant goes to zero for r tend- 
ing to infinity. This is the proper "clean limit" . In the 
opposite limit of strong disorder, the Gilbert damping 
decreases monotonically to zero. This is a manifestation 
of the D'yakonov-Perel effect: in a strongly disordered 
system the instantaneous axis of spin-precession changes 
too rapidly to allow an effective loss of spin orientation. 

Instead of using the equation of motion to calculate the 
Gilbert damping, we could evaluate directly the imagi- 
nary part of the transverse spin response function. The 
two approaches are equivalent and the details of the cal- 
culations of the transverse spin response function are 
summarized in Appendix D. We will see in the following, 
however, that the treatment in terms of the spin-channel 
conductivities allow us to gain useful intuition for the 
inclusion of e — e interaction effects. 



2. Interacting 2DEG 

The role of the Coulomb interaction in spin transport 
was analyzed in detail in Refs. I30lf3ll . but only in the 
absence of SO coupling. In the spin-polarized electron 
liquid, with identical scattering times for up and down 
spins and a = the result is 



3?ecr s ||(0) = cr D 



1 + p 2 JT 
1 + JT 



(60) 



where 7 is the spin-drag coefficient - the rate of momen- 
tum transfer due to the Coulomb interaction between up 
and down spins. The derivation of Eq. (|60[) is presented 
for completeness in Appendix B. 

Notice that the effect of the spin Coulomb drag van- 
ishes (as expected) when the electron gas is fully spin 
polarized, i.e. when p = 1. This is because in this 
limit there are no minority spin carriers to exert a drag 
on the majority spin carriers. More generally, the form 
of Eq. ([BT?|) can be understood in terms of the coupling 
between charge and spin currents in a spin-polarized 
electron gas. A "spin electric field" E s = = —E\ 
drives not only a spin current, but also a charge current 
j c = par>E s (notice that 7 does not enter here). The 
charge-current acts on the spin current as an additional 



spin electric field — JTp 2 E s . The spin current 

then responds in the usual way (i.e. with the unpolar- 
ized spin conductivity u s \ p= q = to an enhanced 

spin electric field E s (l + -/p 2 r). This gives Eq. (|60|) . 

The low-temperature behavior of 7 is approximately 
given by the formula 



7(T) = 



E F 

h 
1 



E F J 6(1 



_ p 2)3/2 

dx 



V 2 (!-P) 



) 2 vt^v / t 



ax 



(61) 



where a = (1 — p)/(l +p).— For zero polarization, 7 di- 
verges logarithmically as can be easily by putting p = 
under integral in Eq. (|6ip . Thus 7 scales quadrati- 
cally with temperature for T/Tf <C 1, but has a non- 
monotonic behavior as a function of polarization. 




FIG. 5: The behavior of renormalization factor (l+p 2 7r) /(1 + 
77") associated with electron-electron interactions as a func- 
tion of polarization for various r s . 

The behavior of the spin conductivity renormalization 
factor (1 +p 2 7-r)/(l + jt) of Eq. ([60]) is shown in Fig. 5 
as a function of p for several values of r s , assuming 



Epr/h = 10/rf and k B T/E F 



! J10. We see that 



the reduction of the spin conductivity and hence of A 
due to Coulomb interaction is significant for small and 
intermediate polarizations, especially at larger values of 

The calculation of the spin drag coefficient in the pres- 
ence of SO coupling poses, of course, a more difficult 
problem. However, in view of the fact that the SO en- 
ergy scale is usually much smaller than the Coulomb 
energy scale such a detailed calculation is not urgently 
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needed. Coulomb interaction corrections are adequately 
taken into account by multiplying the non-interacting re- 
sult of Eq. ([55]) by the correction factor on the right hand 
side of E. (JBDJ. This gives 



A = A. 



1 + p 2 JT 
1 + 7T 



(62) 



Eq. (|62|) has the interesting property that it scales as r 
for r jp 2 T 3> 1 and 77 > 1 as well as in opposite limit 
when 7p 2 r <C 1 and 7T <C 1. The e — e interactions plays 
a role mostly in the regime 1 < 77 < 1/p 2 which shrinks 
to zero as p increases. The behavior of Gilbert constant 
as a function of 7T and polarization p is presented in 
Fig. 6 and Fig. 7. For strong polarization the factor (1 + 
-P 2 7 r )/(l+7 r ) ~ 1 an( i the Gilbert damping scales with a 
scattering time r. The situation looks different for a very 
weak polarization (see solid line in Fig. 7): the Gilbert 
damping scales as t / (1 + 7-7-) for a weak e — e interactions 
but for strong e — e interactions it saturates to a constant 
value proportional to I/7. 

The dependence of A on polarization is more compli- 
cated. This is a consequence of the fact that the spin- 
Coulomb drag depends on polarization. For weak e — e 
interaction A increases linearly with p for small p and sat- 
urates for large p. For strong e — e interactions (7T > 1) 
A increases linearly with p for small p, and increases even 
faster than p for strong polarization. 
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FIG. 6: Gilbert damping constant as a function of 7-r and 
polarizatioon p for nio = 10 12 cm~ 2 and ra e // = 0.067m e . 
We chose the limit of weak spin orbit coupling ak,F /<^o = 0.02 
and T/T F = 0.02. 




10 100 



FIG. 7: The Gilbert damping constant as a function 7-r, the 
two-dimensional cross-section through the 3D plot presented 
in Fig. 5. 



B. Anisotropic case 

Let us now consider the general anisotropic case, in 
which the direction §3 of the magnetization forms an 
angle 6 with the direction e z of the electric field that 
generates the SO effective magnetic field. An exact 
calculation of the torque-torque response function is 
complicated, even in the nonintcracting case, by the 
lack of rotational symmetry about the §3 axis. How- 
ever, in the limit of small a, the spin degrees of free- 
dom are decoupled from the orbital degrees of free- 
dom and we can take advantage of rotational symmetry 
about the 63 axis in spin space to simplify the calcula- 
tion. What happens is that response functions such as 
((J2n Sn3Pnx;J2 n S n ±Pnx)}u> vanish by symmetry when 
S± denotes one of the two transverse components, S± or 
§2, of the spin operator. On the other hand, the response 
function ((J2 n S n ±p nx ; Y, n S n ±p nx )) u does not depend 
on which component S n j_ we choose to consider. 

Going back to Eqs. (|32II33[) and making use of Eq. (|4"Tj) 
we see that, to second order in a we have 



A, 



Hp 



3fo7 B ||(0) ^2 e fte<T a x(0) 



CD 



0\D 



(64) 



It is also possible to include the frequency dependence 
of the spin Coulomb drag correction through the replace- 
ment of i by i — iia. A straightforward calculation gives 



A(w) = \ ni {bj) 



(1 + 7t)(1 +p 2 7-r) + {lot) 2 
(1+7t) 2 + (wt) 2 



(63) 



and 
A 2 



a 2 E F r 

hp 



cos O h sm O 



0\D 



(65) 

where the transverse spin- channel conductivity o~ s ± is de- 
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fined, in analogy to Eq. (|44j) as 



V. GENERALIZATIONS AND CONCLUSION 



nx u 

Wecr s ±(uj) = —~sm- 



m 2 V 



h 2 i 



(66) 

The calculation of a s ± in the absence of SO coupling 
and e — e interactions is straightforward and yields 



3te<r s± (0) 



1 + (u Tf 



(67) 



Notice that this differs from the longitudinal spin-channel 
conductivity cr s || , calculated in the same approximation, 
simply by the factor [1 + (luot) 2 ]" 1 . This takes into ac- 
count the non-conservation (precession) of the transverse 
spin caused by the magnetic held in the §3 direction. 
Therefore, in the limit of weak SO coupling and no in- 
teraction we get 



Ai 



of Eft 
hp 



1 



sin 2 



1 + (loot)' 



and 



a 2 E F r 

hp 



cos 2 



sin 2 

1 + (uoTf 



(68) 



(69) 



The effect of the spin Coulomb drag can be included, 
at the same level of approximation, by multiplying the 
two spin-channel conductivities by appropriate renormal- 
ization factors. Namely, we use Eq. (|60|) for cr s ||(0), and 



$tea s ± (0) = (jd 



1 



■ 1±T 



(loot) 2 + (1 + 7± r) 2 



(70) 



where 7_l is a transverse-spin analogue of the longitudinal 
spin Coulomb drag coefficient, as discussed in Appendix 
C. 

An unsatisfactory feature of Eqs. (f6"B")) and (|6"9"|) is that 
they are valid only for akp <C ujq. Notice however that 
the angular dependence of the Gilbert constants in these 
equations reduces to simple angular factors multiplying 
the isotropic result. This suggests that we take care of 
the problem simply by multiplying the full isotropic non- 
interacting result for A(w), given by Eq. ([55]) . by the same 
angular factor. The resulting expressions: 



Ai — X r 



5Recr s ||(0) 



(TD 



and 
A 2 



A, 



2 ^^11(0) 

cos 

(TD 



2 Q aic(7 s j_(0) 



sin 



<TD 



(71) 



(72) 



with A„i, 9?eo- s ||(0) and 3?ecr s _L (0) given by Eqs. (55]), (|60j) 
and (|70[) respectively, vanish for loq tending to zero, and 
should be reasonable not only for akp *C but also 
for ak,F 3> luq provide both energies are small compared 
to the Fermi energy. Furthermore these formulas can 
be generalized to finite frequency simply replacing A„i 
by A„i(w) (Eq. ([59")) ) and replacing - by - — ilo in the 
Coulomb rcnormalization factors. 



In this paper we have shown that (1) the Gilbert damp- 
ing constant for a homogenously magnetized electron gas 
with SO interactions can be exactly expressed in terms 
of the torque-torque correlator which comes from SO in- 
teractions; (2) In the special case of the Rashba model 
the Gilbert damping can be expressed in terms of the 
spin-channel conductivity. Based on this connection, we 
have discussed the behavior of the Gilbert damping con- 
stant in a two-dimensional electron gas with Rashba SO 
coupling as a function of magnetic held, SO interaction, 
e — e interaction, and disorder. These calculations, while 
based on linear response theory, do nevertheless provide 
the input for the nonlinear LLG equation, provided it is 

recognized that A depends on the instantaneous orienta- 
tion of the magnetization. 

It should be clear that point (1) above is completely 
general, while point (2) depends on a specific feature of 
the Rashba spin-orbit interaction, namely, the fact that 
the electric field is independent of position. This poses 
the question: to what extent is the connection between 
Gilbert damping and spin conductivity transferable to 
general, non-homogeneous, spin-orbit coupled systems? 

In this conclusion we argue that the connection is in- 
deed broadly applicable to itinerant-electron ferromag- 
nets, under assumptions similar to those that are used to 
justify the local density approximation (LDA) in the den- 
sity functional theory of electronic structure. We recall 
that in a multi-band all-electrons theory the SO interac- 
tion has the form: 



ks ° = ~ 2^ £ ( W ( r ") X P") ' S» 



(73) 



where V(r) is the self-consistent Kohn-Sham potential, 
typically given by the local density approximation. The 
gradient of V(r) defines a privileged direction at each 
point in space, which we denote e z (r). Then the spin- 
orbit interaction can be recast in the form 



H 



so 



E¥[^») x P»]-S», (74) 



where a(r) is a position dependent SO coupling constant, 
which reflects the local magnitude of the electric field 
seen by the electron, and e z (r) is the local direction of 
this field. Another privileged direction is the direction of 
the local magnetization, which we call §3 (r). Thus we are 
back to our model hamiltonian (f55|) (we already remarked 
that the two-dimensionality of the model is not essential 
to our treatment) with the crucial difference that a, e z , 
and §3 are functions of position^ 

Now in the spirit of the LDA let us assume that the 
Kohn-Sham potential, its gradient, and the direction of 
the magnetization are all slowly varying on the appropri- 
ate microscopic length scale (e.g., the electron-electron 
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distance, or the size of the unit cell, or the scattering 
mean free path). Then we can equiparate each small 
volume clement of the system to a homogeneously mag- 
netized electron gas of density n and polarization p point- 
ing in the 63 direction, with a SO coupling constant a(r) 
and a local anisotropy axis e z . With this identification, 
all the results obtained in the previous section become 
immediately applicable to each volume clement. 

An obvious objection to this procedure is that the elec- 
tronic density of real materials is not slowly varying on 
the atomic scale. In spite of this, however, it is well 
known that the LDA works very well in materials be- 
cause the exchange correlation potential is controlled by a 
spherical average of the exchange-correlation hole, which 
is reasonably close to the spherical hole of the homoge- 
neous electron gas. We believe that a similar averaging 
may also work for the magnetization dynamics. Then the 
Gilbert damping would be given by the formulas derived 
in this paper, only with the appropriate coarse-grained 
values of a, e z , and §3. 

VI. ACKNOWLEDGEMENTS 

We thank Allan H. MacDonald and Arnc Brataas for 
enlightening discussions. This work was supported by 



NSF Grant No. DMR-0313681 and supported in part 
by the National Science Foundation under Grant No. 
PHY99-07949. 



VII. APPENDIX A: LONGITUDINAL 
SPIN-CHANNEL CONDUCTIVITY FOR 
NONINTERACTING ELECTRONS 



In the calculations that follow we set h = 1 for 
convenience. To obtain the real part of longitu- 
dinal spin-channel conductivity we need to calculate 

SB! " Sn3 P nX '^p™ Snafinx))^ ^ a ma g ne ti za ti on p ar _ 

allel to the z axis (i.e. 6 = 0) the above expression 
simplifies to Sm (( SzP *'Sz p *))^ wne re S z is the z compo- 
nent of the spin operator and P x is x component of the 
total momentum. For non-zero uj we have to evaluate the 
following integral: 



where 



G(p, uj) = — (75) 

uj - e p — h p S + i/2Tu 

is the disorder-averaged Green's function near the Fermi level, e p is the kinetic energy relative to the Fermi level, 
hp = (ap y , —ap x ,u)o) is the effective magnetic field, S = ct/2, where a is the vector of Pauli matrices and 

— = -signH . (76) 



Using the fact that the integration over energy involves 
only the states around the Fermi energy, i.e. e + ui > Ep 
and e < Ep allows us to integrate over e and cancel out 
the 1/ui on the r.h.s. of Eq. (|75|) . Then one can see 
that the formula for a non-zero ui can be obtained from 
the uj = limit with the substituition 1/r — > 1/r — iuj, 
leading to Eq. (jST)]) . 

In the ui = limit G(e + uj) and G(uj) should be sub- 
stituted by G(0 + ) and G(0~), respectively. The ladder 
vertex corrections are found by solving the self-consistent 
integral equation for A. For an electron-impurity poten- 



tial of the form U 5(r) this equation is 

K = pS z + U 2 J 0^G( P ', 0+)AG(p', 0") (77) 

Its solution gives the vertex correction of the following 
form: 

where S + = S x + iS y and U 2 = 1/mr. 

The calculation of the spin-channel conductivity in 
Eq. (|75"|) without the vertex corrections (a single bubble) 
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gives: 



((S Z P X ;S Z P X )), 



vnnT ( n 2 T 2 cosS + 1 

n 2 r 2 + 1 



(79) 



Substituting the elements of the resistivity matrix into 
Eq. (|54"|) one obtains: 

[lo 2 + l/r 2 (l + 7rp 2 )](7 + 1/t) - 1P 2 lo 2 



while the final result including the ladder vertex correc- 
tions is: 



((S Z P X ; S Z P X )) U -mnr cos 2 S[l + ±ft 2 r 2 (l + cos 2 5)] 
wV ~ 4 cos 2 5 + \VL 2 t 2 (I + cos 2 5) 2 



Using Eq. and the definition of spin channel con- 

ductivity we finally obtain the formula Eq. (|5"Tj) for the 
real part of longitudinal spin-channel conductivity: 



5ffei7 s || (0) = (Jd cos 2 5 



l + i(ftT) 2 (l + cos 2 5) 
cos 2 <5+ i(ftr) 2 (l + cos 2 8) 2 



The single bubble calculation recovers the behavior of 
the longitudinal spin-conductivity for weak SO interac- 
tions: Eqs. ([79")) and (|5D|) coincide for 5 = 0, i.e. fora 
zero SO coupling. However, for strong spin-orbit inter- 
actions the vertex corrections are absolutely necessary. 
Moreover, only Eq. (|8D|) gives the correct result for zero 
magnetic field, i.e. 5Recr s i[ (0) = 0. 



VIII. APPENDIX B: LONGITUDINAL 
SPIN-CHANNEL CONDUCTIVITY IN THE 
PRESENCE OF e — e INTERACTIONS 

In the absence of SO interactions the structure of the 
2x2 resistivity matrix can be deduced from the equation 
of motion for electrons with spin a and velocity v a : 



-eN„E a + Fg 



(81) 



where m* is the effective mass of the electrons, a = —a, 
F^ s is the average Coulomb force between electrons of 
opposite spin orientations, which is proportional to the 
difference of their drift velocities. The Fourier transfor- 
mation of Eq. (f8"Tj) leads to the following formula for the 
spin current: 



a s u (oj) = <J D 



(lo 2 t+1/t)[lu 2 + (1/t + 7 Y 



If ui = 0, Eq. |Jg5J) simplifies to Eq. (|50J> : 

1 +p 2 ^T 



(85) 



3?ecr s ii (0) = ctd- 



1 + 7T 



For zero polarization the spin and charge channels are 
decoupled and the inverse of the effective scattering time 
consists of two contributions, one connected with disor- 
der and second one associated with e-e interactions: 



(86) 



For nonzero polarization the spin and charge channels 
are mixed and the additional term 7p 2 r appears in the 
numerator. 



IX. APPENDIX C: TRANSVERSE 
SPIN-CHANNEL CONDUCTIVITY IN THE 
PRESENCE OF e e INTERACTIONS 

In this appendix we use the semi-classical approach 
(similar to the one presented in a previous appendix) 
to calculate the real part of the transverse spin-channel 
conductivity. According to Eq. (|66|) : 



^iea xx = $lea± = 



mrV 



(82) 



so we are looking for a spin current response on the 
perturbation p x S x . This perturbation is generated by 
the SU(2) vector potential A X S X where S x is the x- 
component of the spin operator. As a consequence one 
can show that the "electric field" that drives the change 
of spin current does not depend on spin (see Eq. flSTJ be- 
low). We neglect SO interactions in this calculation and 
work in the circularly polarized basis in which the spin- 
conductivity is diagonal and has two components af{uj) 
and aj(cu). The diagonal character of the conductiv- 
ity in this basis is equivalent to the statement that the 
current of right-handed circularly polarized electrons, j + 
does not interact with the current of left-handed circu- 
larly polarized electrons j_. The semiclassical equation 
of motion for these currents is 



Using Eq. (|82)) . we can show that the resistivity tensor 
has the form: 



Tf n 1 
— -^C'V 



H? 1 V 



^t7 

-iuj + J- + m 7 



eS i± v ix -/± 
(87) 



Inverting the resistivity matrix we get the following for- 
mula for the spin-channel conductivity a s \\ : 



where S± = S x ± iS y , v x is the a;-component of the ve- 
locity and 7j_ is in plane spin-Coulomb drag coefficient. 
The Fourier transformation of Eq. (|87|) leads to the fol- 
lowing formula for the left and right circularly polarized 
spin current densities: 



o" s ||H = 



AU +PTT + PU +P1T 



(84) 



-ne 



E 

m* 



,3±,. 
H 1- J±7-L 



(88) 
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where we used the fact that in the absence of SO coupling: 
dS± 



dt 



= —i{S±,H] = ±iu>oS± 



(89) 



Then the circularly polarized spin conductivities have the 
form: 



ne 



1 



(90) 



where we assumed that the spin-Coulomb drag coefficient 
is the same for both circular polarizations^ Accordingly, 
the real parts of spin conductivities have the following 
form: 



3?e(4(w) 



7_l + 1/t 



m* (wtw ) 2 + (7i + 1/t) : 



(91) 



Assuming that u» 1/r, 7, u>q we find the real part of the 
transverse spin channel conductivity <j^_(uj) as follows: 



Eq. (EL 



5fte<4H + $eaL(uj) _ ne 2 7± + l/i 



2to* 



.(92) 



Eq. (I92[) separates the disorder and interaction terms. 
Since we have already included disorder in the phe- 
nomenological Eq. (|88p . it is justified to derive the trans- 
verse spin-drag coefficient 7j_ by comparing Eq. (|92p 
with the transverse spin-channel conductivity found from 
Eq. (f66|) . Using two times the identity (fl"9|) we can rewrite 
Eq. (|66p in terms of the force-force correlation function: 



Rea±(uj) = 2T7^ m 



m 2 V 



(93) 

Then by comparing Eq. (|92[) with Eq. (|93[) we obtain in 
the clean limit: 



7_l(w) = -- 



nm*V u> 
f°° ^ qq 3 h ,2 erPuJ ~ 1 /"" d ^,Q ! "ix„„(q,^') Sm Xx K ('7,w - w') 



2-K 2 nm* Jo w Jo (e-/ 3 "' - l)( e -/3(^-^') - 1) 

kp f°° 3 2 /"" e ~ 2l3EFU: - 1 j _ / 3mx„™(g, u)')^smxxx(q, to - ui') 



d ™ v « I r, ' (94) 



2-K 2 nm* J q Jo u ( e -2/3E P w' _ ^^-2^(0-0') _ x) 

r 



where q = q/kF, to = Hlj/2Ef- For small 

V 9" 4 

and 3mx„„(g,u) = 3mxt(g,w) + 3roxx (g, w), where 
yq^m X ] {q, w) = — /= — -= (96) 



and Smxi (<7> <^0 i s obtained from the above equation sim- 
ply by changing the sign of p. Vq is the dimension- 
less screened interaction potential and r s = — is 
the dimensionless Wigner-Seitz radius. In the limit of 
T/Tp <C 1 this formula can be evaluated analytically 
leading to the following integral on the transverse spin 
Coulomb drag coefficient: 



where 



7J-(0) = ^ (^) 2 [F{r.,p) + F(r s , -p)} (97) 

f b{p) qdq 1 
F(r s ,P) = / ,_ q l , 2 7= — , (98) 

a( p) (g + \/2r s )^ 



for 6(p) > a(p) and zero otherwise, and the limits of integration are 



o(p) = IVl+p- VI -p| , 6(p) = min(VT + V + \/l - P, 2\/l + p) . (99) 

I 



For zero polarization the transverse spin-Coulomb drag the isotropy of the electron gas. The behavior of the 
coefficient coincides with the longitudinal one because of 
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FIG. 8: The behavior of ^y±r as a function of polarization for 
various r s . We assume Eft/H = 10/r2 and ksT/EF — rf/10 



transverse spin- Coulomb drag coefficient is presented in 
Fig. 8 for various values of r s . The effect of the electron- 
electron interactions is stronger for larger r s and there- 
fore 7j_t is larger. For a given r s , "f± decreases with 
increasing polarization and saturates for large polariza- 
tions to the value of jt± corresponding to the value of the 
F(x,p) integral for p = 1 (see Eq. (|98|) ). The downward 
step at polarization p — 0.8 is due to the fact that the 
limits of integration for F(r s , —p) come together at this 
value of p and furthermore the integrand has a singularity 
that causes the integral to drop to zero discontinuously 
at this point. 



X. APPENDIX D: DIRECT MATSUBARA 
CALCULATION OF THE NONINTERACTING 
TRANSVERSE SPIN RESPONSE 

The Gilbert damping constant is directly related to 
the imaginary part of the transverse spin susceptibility, 
which docs not depend on the equilibrium parameters 
(for example magnetization) of the system. It is useful 
to show that the calculation of the transverse spin sus- 
ceptibility via the torque-torque correlator and the more 
direct calculation by Matsubara approach arc equivalent. 
We consider the same Hamiltonian as in Eq. ([55)) . The 
transverse spin response x h is defined as follows: 



X-+(r,r';t,t') 



2V 



Q(t-t')([s_(r,t),s + (r',t'))) 



(100) 

where s(r) = Tr[(T/5(r)/4] is the spin density (s± — s x ± 
is y ) in terms of the density matrix p Q ^(r) = ^'j 3 (r)\l/ Q ,(r). 



It is straightforward to show that the disorder-averaged 
retarded single-particle Green's function in the represen- 
tation of Rashba subbands is given by 

^(k,k';.) = 5k , k ,L> k ( G ^ k ^) c^)^ 

(101) 

(102) 



where 



G*(k,u,) = 



lo - (e k ± Cl) + i/2r ' 



e k = fr 2 k 2 /2m — pi , and t/ k are spin-rotation matrices: 



where 
and 



cos(,5/2) isin(<5/2) e - 8¥ "< 
-isin(<5/2)e i¥,k -cos(<5/2) 



w 



^ y/p 2 + a 2 



n = 



(103) 



(104) 



(105) 



arc defined as in main text. Note that, in general, G A = 
(G R y , in order to obtain the advanced Green's function. 

We do not assume any a priori hierarchy in the three 
energy scales H, akp, and 1/r, but we consider for sim- 
plicity the limit e/ej? 1, where e is any of the men- 
tioned energy scales. Disregarding terms of order e/ep 1 
k can be set to &f m the relevant expressions. 

It is most convenient to calculate the spin response 
function in the Matsubara formalism. For a uniform per- 
turbation, one gets in general 

X -+{iVl n ) = -— G TT (k',k;iw m ) 

k,k';m 

x Gn(\s.,\s!;iw m +iO„) , 

(106) 

where uj m = 27rT(m + l/2) and Cl n = 2-nTn are fermionic 
and bosonic Matsubara frequencies, respectively. T is the 
absolute temperature and V the total volume of the sys- 
tem. The desired retarded response function is obtained 
by analytic continuation ift n — > to + i0 + . 

It is important to take into account vertex corrections 
when averaging Eq. (|106[) over disorder. Otherwise, the 
response would have a low-frequency dissipativc compo- 
nent even in the absence of SO interaction, i.e., a = 0, 
if the disorder-averaged Green's functions where inserted 
in Eq. (|106[) . Namely, we find without vertex corrections: 



g 2 m*T 

v n^r 



i 



i 



sin 2 S 



(fir) 2 

(107) 

which diverges in the clean limit, r — > oo, even when 
there is no SO interaction, i.e., 5 — 0, and vanishes in 
the dirty limit, r — > 0. 
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Taking into account vertex corrections is conceptually 
straightforward but technically somewhat tedious. Defin- 
ing the tensor II: 



n a! ' / 3 iCt/ 3'(iO„) = ^ G a >p(k!,k.;i(j rn ) 



kk' m 



x G a p>(k,k';iu; m + in n ), (108) 

whose TT) II component determines the response 

function according to Eq. (|106|) . namely, x h = 

— (T/2T^)II|| i |^,we find by summing the ladder diagrams 

n a '0 !a 0'(£O„) = ^r Q ,/ ij( 3'(i^ n ,iw rn ) 
m 

x T^ ;QJ (ifi n , iwm) , (109) 

where 

TC a '/3, Q /3'(«r2ri, iw m ) = J^Gq.'^ (k; iw m ) 

k 

xG a/ 3'(k;iw m + iO n ), (110) 

is the uj m contribution to the II tensor without the vertex 
corrections and 

1 



2lTVT 



x Ti/3 !a j(iQ n ,iu m )T a ' it jff'(ifl n ,iu} m ) (111) 

is the recursive relation for the vertex tensor T. v = 
m* /2tt is the density of states per spin. If a tensor A 
can be decomposed into matrices B and C by A a ip^ a p> = 

Ba'pCafS', I will write A = A ® B. All subscripts in 
Eqs. (fill))) . (fTTU]) . and (TlTT]) are spin-1/2 indices taking 
values | and J. , and the summation over repeated indices 
i and j is implied. 

In order to calculate the II tensor (|109[) . we first fol- 
low the standard procedure of replacing the sum over 
m by the integral J c dz ta,nh(z / 2T) with the contour C 
around the imaginary axis. The contour is then deformed 
towards the ±oo along the real axis, around the branch 
cuts at Imz = and iCl n . Performing integrals over the 
branch cuts, we get 

nM = £^tanh(^) P^(e,e + W ) 
-X A ' R (e, e + u) + x A ' R (e - u, e) - x A ' A (e - u, e) 



(112) 



where 



X (e,0 = Xo ( £ > e ) 



£]T G A "(k, e)x X ' y (e, e')G r (k, e') (113) 



is the recursion relation corresponding to Eq. (jllip 
(X, Y = A or i? and = 2twtV0 and 

X^ Y (e, e') = G A >, e ) ® G Y (k, e') . (114) 



Since we are interested in dissipation, let us define 
2ImTI = -i(n-nt), where {W} a ^^ = H* 0a , >0 , a . We 

then take the following steps in order to evaluate Imfl: 
First, Eq. (|113[) is iteratively expanded in terms of the 
disorder-averaged Green's functions G x , then the "un- 
rolled" x X ' Y is substituted into Eq. (pJJ) , Differen- 
tiating the resulting expression with respect to to and 
taking the imaginary part, the integral can be trans- 
formed integrating by parts into J c?etanh(e/2T)3 e ... = 
— J ded e tanh(e/2T)..., and, assuming low temperatures, 
we approximate d t tanh(e/2T) w 25(e), so that the dis- 
sipation is naturally governed by electron-hole pair exci- 
tations near the Fermi surface. The infinite summation 
series is finally "rolled" back into a recursive relation and 
we obtain the following expression for the spin response: 



d u Imx-+(u))\ u=0 = - 



g 2 m* 



V 8tt 



*R,R 



*A,R 



X 



R.A 



'TT,U 



(115) 



where all %'s are now evaluated at e, e' = 0. 

The problem is thus reduced to calculating x X ' Y (0, 0) 
using Eqs. (| 1 13|) and (|114[) . It is still somewhat tedious 
but now totally straightforward. We first do the angular 
integration fixing the absolute value of k in Eq. (|114[) . 
This averages over the rotation matrices ()103|) . The re- 
sulting angle-averaged tensor (|114[) then decomposes into 
6 components: 1 ® 1, 1 ® g z , g z ® 1, g z ® g z , g + <8 
and a- ® g + (the last two arc solely due to the SO cou- 
pling), where g± = g x ± io"^. The respective prefactors 
are given by integrals over the absolute value of momen- 
tum near the Fermi energy, which are trivial to evaluate 
after linearizing the dispersion at the Fermi level. We 
then make an ansatz that the tensor x > (0,0) deter- 
mined by Eq. (fTT3|) can also be expanded in terms of 
the same 6 components and, after plugging the expanded 



X 



X ' Y (with unknown coefficients) and Xo' Y (with calcu- 



lated coefficients) into Eq. (| 1 1 3|) . we obtain a linear sys- 
tem of equations for the unknown coefficients that deter- 
mine X X ■ Solving this (with Mathematica), we find a 
solution (validating the ansatz), and the dissipative part 
of the spin response is finally given by 



d UJ Imx-+(uj) 



g 2 m*T 



r 2 5 



1 + ±(»t) 2 (1+cos 2 5) 

X cos 2 5+i(Or) 2 (l+cos 2 (5) 2 ' 

(116) 

Using Eq. (|2"6"|) and Eq. ([3H)) one can show that Gilbert 
damping has a following form: 

X = gM V Hm gffiztM - 



w^o xl(u = 0) 
E f t sin 2 5 l + i(fV) 2 (l + cos 2 <5) 
P H cos 2 (5+i(fiT) 2 (l + cos 2 5) 2 

in agreement with Eq. (| 5 5 1) . 



(117) 
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